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Abstract. We report on the hrst multi-epoch, phase referenced VLBI observations of the H2O maser emission 
in a high-mass protostar associated with a disk-jet system. The source under study, IRAS 20126+4104, has been 
extensively investigated in a large variety of tracers, including H2O maser VLB A data acquired by us three 
years before the present observations. The new findings fully confirm the interpretation proposed in our previous 
study, namely that the maser spots are expanding from a common origin coincident with the protostar. We also 
demonstrate that the observed 3-D velocities of the maser spots can be fitted with a model assuming that the 
spots are moving along the surface of a conical jet, with speed increasing for increasing distance from the cone 
vertex. We also present the results of single-dish monitoring of the H2O maser spectra in IRAS 20126+4104. 
These reveal that the peak velocity of some maser lines decreases linearly with time. We speculate that such a 
deceleration could be due to braking of the shocks from which the maser emission originates, due to mass loading 
at the shock front or dissipation of the shock energy. 
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1. Introduction 

Several molecular species are known to exhibit maser emis- 
sion from the sites of newly born stars (see e.g. Elitzur 
1992). This applies especially to the regions where high- 
mass stars (i.e. stars in excess of ~8 Mq) form. In par- 
ticular, the water line at 22.2 GHz, originating from the 
616 — * 523 transition, is the most powerful emitter, attain- 
ing flux densities up to 10 6 Jy in a ~1 km s _1 wide line. 
Since the emitting region (named "spot" ) of a maser line is 
as small as 1 mas or less, H2O masers are excellent tools for 
high angular resolution and low sensitivity observations 
such as those performed with very long baseline interfer- 
ometry (VLBI). These characteristics make it possible to 
use H2O masers as "test particles" to trace the velocity 
field of the gas in the densest and hence most obscured 
portions of molecular clouds. This is of great help in stud- 
ies of massive star forming regions, which are hindered by 
two problems: massive stars are much more distant (a few 
kpc) than low-mass stars and are born deeply embedded 
in their parental cores. Both shortcomings are overcome 
by H2O maser observations: on the one hand, VLBI ob- 
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servations can attain angular resolutions of ~1 mas, cor- 
responding to ~1 au at 1 kpc; on the other hand, obser- 
vations at 22 GHz may easily penetrate even the highest 
column densities. 

For these reasons, H 2 masers have been extensively 
used to pin-point the location of young stellar objects 
(YSOs) and investigate the structure and kinematics of 
the gas in their surroundings. Although some authors have 
suggested that such masers might originate in circumstel- 
lar disks (see e.g. Torrelles et al. fTMR [TQM Goddi et al. 
2004} , the common belief is that they are associated with 
outflows (Felli et al. I1992f) . In one case (Torrelles et al. 
2001J, the maser spots have been found to describe a per- 
fect circle (to one part in a thousand), suggesting that they 
could originate at the interface between a spherical stellar 
wind from a deeply embedded YSO and the surrounding 
gas. 

With this in mind, Moscadelli et al. (2000; hereafter 
MCR) performed 1 mas resolution observations with the 
very large baseline array (VLB A) of the H2O maser emis- 
sion in the well-known massive YSO IRAS 20126+4104, 
which is the best example to date of a high-mass pro- 
tostar associated with a Keplerian disk and bipolar out- 
flow/jet (Cesaroni et al. 113371 IT3331 12TJU51 Zhang et al. 
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[TUMI Hofner et al. H555l Shepherd et al. 12009 • The scope 
of the observations was to establish if the masers were 
co-rotating with the Keplerian disk or expanding with 
the outflow /jet. The former possibility was ruled out by 
MCR, who demonstrated that the H2O maser spots are 
distributed along the axis of the outflow /jet and may be 
adequately fitted with a conical model implying an expan- 
sion velocity of ~23 km s _1 . Such findings are of great 
importance for a better understanding of the structure of 
the jet in IRAS 20126+4104, which is known to be under- 
going precession (Shepherd et al. EOOO). Consequently, to 
establish the current direction of the jet it is necessary 
to observe it as close as possible to the powering source 
(see Cesaroni et al. 1200 5 for a thorough discussion of this 
topic) , which may be best achieved with a VLBI study of 
H2O masers. 

Notwithstanding the results obtained by MCR, their 
study was limited to the l.o.s. component of the velocity 
of the spots. This is a serious limitation in the case of 
IRAS 20126+4104 because the jet axis lies very close to 
the plane of the sky (Cesaroni et al. I1999|) . thus making 
an accurate estimate of the expansion velocity very diffi- 
cult. A small error on the estimate of the inclination of 
the velocity vector may be reflected as a large uncertainty 
on the true expansion speed. In order to overcome this 
problem and hence obtain a more accurate description of 
the jet, we have performed new multi-epoch observations 
with global VLBI, which resulted in proper motion mea- 
surements of the H 2 maser spots. In the following we 
present the results of this study and improve the model 
fit proposed by MCR to take into account the spots' mo- 
tion in the plane of the sky. Finally, we present the results 
of single-dish monitoring of the H2O maser emission from 
IRAS 20126+4104 over a period of ~15 years and interpret 
these in the light of the VLBI findings. 

2. Observations and data reduction 

IRAS 20126+4104 was observed using the Global array 
(including a subset of EVN 1 antennae plus the 10 an- 
tennae of the VLBA 2 ) for 18 hours at three epochs, on 
November 9 and 26, 2000, and on March 1, 2001. A 
group of 6 EVN antennae (Effelsberg, Jodrell, Medicina, 
Noto, Onsala, and Shanghai) took part in the observa- 
tions at each epoch, whereas Metsahovi and Robledo ob- 
served only the first two epochs and the second epoch, 
respectively. The observations were performed in phase- 
referencing mode alternating scans on the maser target 
with scans of the phase-reference source J2007+4029, with 
a cycle time of ~1 min. J2007+4029 is an intense quasar 
belonging to the ICRF (International Celestial Reference 
Frame) list, with a relative separation from the target 

1 The European VLBI Network is a joint facility of 
European, Chinese, South African and other radio astronomy 
institutes funded by their national research councils. 

2 The VLBA is a facility of the NRAO, which is operated by 
Associated Universities, Inc., under contract with the NSF. 



maser source of 1?5. For the purpose of bandpass and 
phase calibration, a scan of a few minutes on one of the 
strong, compact calibrators J2007+777 and J2002+4725, 
was observed every 1-2 hours. 

The stations recorded an aggregate of 16 MHz band- 
width in both circular polarizations, centered at the LSR 
velocity of —3.5 km s -1 (based upon a rest frequency 
of 22235.0798 MHz), using 2-bit sampling. The correla- 
tion was made at the VLBA correlator in Socorro (New 
Mexico) using 1024 spectral channels uniformly weighted, 
which gives a channel separation of 0.21 km s _1 . 

The amplitude and phase calibration of the observed 
visibilities and the mapping were accomplished with the 
NRAO AIPS package, following the same procedure de- 
scribed by MCR. At each epoch, the maser phase-reference 
channel was chosen to be the one with the strongest emis- 
sion, whose velocity changes slightly across the epochs 
(from —14.7 km s _1 to —15.5 km s -1 ), and whose struc- 
ture consists of a single, unresolved spot. The absolute po- 
sition of the reference maser spot has been derived from 
the maps produced after applying the corrections calcu- 
lated using the phase-reference source data. These cor- 
rections were obtained in two steps. First, after applying 
the calibrator corrections, the visibilities of the phase- 
reference source were fringe fitted to find the residual 
fringe rate produced both by differences in atmospheric 
fluctuations between the calibrators and the source, and 
by errors in the model used at the correlator. Then, af- 
ter correcting for the residual fringe rate, the visibili- 
ties of the phase-reference source were self-calibrated to 
remove any possible effect induced by extended spatial 
structure. For the first observing epoch, the absolute po- 
sition of the reference maser spot is: RA(J2000) = 20 h 
14 m 26 s .0253, Dec(J2000)= 41° 13' 32'/666, estimated to 
be accurate within a few tenths of milliarcsecond. As a 
check, we have also calculated the absolute position with 
the reverse phase-referencing technique (i.e. mapping the 
phase-reference source after applying the corrections eval- 
uated working with the maser data), and found fully con- 
sistent results. 

Using the AIPS task "IMAGR" , l"xl" (E x N) ta- 
pered, channel-averaged maps centered on the reference 
spot were produced for the velocity range from —28 to 
+3 km s _1 , which includes all the emission apparent in 
the total power spectra. We note that, at each epoch, 
the maser emission extends over a velocity range nar- 
rower than that of the VLBA observations by MCR, for 
which spectral features in the velocity interval from +3 to 
+15 km s _1 were also observed. The field of view covered 
by the tapered maps is about twice the sky area over which 
water maser emission is detected in the MCR's maps. 
The bulk of maser emission is found within a distance 
of «200 mas from the reference spot, with the exception 
of a single feature detected at a distance of ^530 mas. 
In order to map the maser emission at full angular and 
velocity resolution, at each epoch and for each velocity 
channel in the range from —32 to +10 km s^ 1 , we pro- 
duced 074 x Of.' 4 (E x N) naturally weighted maps cen- 
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tered on the reference spot. In addition to that, smaller 
maps were also created at the position and velocity of the 
more detached feature. The CLEAN beam was an ellipti- 
cal Gaussian with a full width at half maximum (FWHM) 
size, slightly varying from epoch to epoch, of 0.6-0.8 mas 
along the major axis and 0.5-0.6 mas along the minor 
axis. In each observing epoch, the la RMS noise level on 
the channel maps is close to the theoretical thermal value, 
3 mjy beam -1 , for channels where no signal is detected, 
and increases to 13 mjy beam -1 for channels with the 
strongest components. 

Every channel map was searched for emission above a 
conservative detection threshold (in the range 5-15 a), 
and the detected maser spots were fitted with two- 
dimensional elliptical Gaussians, determining position, 
flux density, and FWHM size of the emission. Hereafter, 
we use the term "spot" to indicate the mean position of 
a collection of spectrally and spatially contiguous maser 
spots. As the typical line width of the water masers is 
greater than 0.5 km s -1 , a maser feature is considered 
real if it is detected in at least three contiguous channels 
(0.21 km s -1 wide) at the same position, within an uncer- 
tainty equal to the FWHM obtained with the Gaussian 
fit. 

The uncertainty in the relative positions of the maser 
spots is estimated using the expression (Reid et al. Q988 ) 

^=jj FWHM (1) 

where FWHM is the un-deconvolved spot diameter, / is 
the peak intensity and a is the RMS of the map evaluated 
over a region where no signal is present. Depending on the 
spot intensity, the relative positional uncertainty varies in 
the range ~ 0.1-100 //as. 

3. Results 

3.1. Maser positions and velocities 

Tabled lists the parameters of the detected maser spots. 
Column 1 gives the spot label number. Columns 2 and 3 
report the positional (RA and Dec) offsets (measured on 
the first epoch of detection) calculated with respect to the 
reference spot (labeled # 1). Such offsets are estimated 
from the (error- weighted) mean positions of the contribut- 
ing maser spots. The positional uncertainties are evaluated 
by taking the weighted standard deviation of the spot po- 
sitions. Columns 4 and 5 list respectively the integrated 
flux density, S v , and the line-of-sight velocity, Vlsr, of 
the spot highest-intensity channel, both averaged over the 
observational epochs for the time-persistent spots. 

The absolute position of a spot at each observing 
epoch is calculated adding the (RA and Dec) offsets to 
the absolute position of the reference spot. The absolute 
proper motions have been calculated performing a (error- 
weighted) linear least-squares fit of the absolute positions 
with time. Columns 6 and 7 of TableHreport the projected 
components, respectively along the RA and Dec axis, of 



the measured absolute proper motions, together with the 
formal errors of the linear least-squares fit. Using a dis- 
tance to IRAS 20126+4104 of 1.7 kpc, Columns 8 and 9 
report the absolute proper motion components after cor- 
rection for parallax, motion of the Sun with respect to the 
LSR, and galactic rotation. 

3.2. The distance to IRAS 20126+4104 

The absolute proper motions measured by us are the sum 
of the intrinsic proper motion due to the velocity of the 
maser spots with respect to the YSO and the apparent mo- 
tion due to the annual parallax and motion of the source 
with respect to the Sun. Therefore, one of the possible 
outcomes of a proper motion study is the derivation of 
the distance to the source by fitting the annual paral- 
lax. However this technique cannot be applied in our case, 
as the few measurements (3) and the short intervals be- 
tween them (<95 days) are insufficient to properly sample 
the apparent position of the spots as a function of time - 
at least 5 points evenly distributed over 1 year would be 
needed (see e.g. Vlemmings et al. 2003). Also, for the H2O 
masers in IRAS 20126+4104 the situation is complicated 
by the fact that the intrinsic proper motion of each spot 
(4-14 masyr -1 ) is greater (see Sect. l3.3|l than that due to 
the annual parallax (^1.2 masyr -1 ). 

The same considerations hold when deriving the in- 
trinsic proper motion of each spot: the apparent proper 
motion due to the parallax and source velocity must be 
subtracted 3 from the observed proper motion. Such an ap- 
parent motion depends on the distance and consequently 
also the estimate of the intrinsic motion depends on it. 
One may thus use this dependency to constrain the dis- 
tance. Two approaches are possible. 

The first uses the fact that, to a first order approxima- 
tion, the maser spots should move along a given direction, 
since they are participating in the expansion of the jet 
from IRAS 20126+4104. This implies that for each spot 
the positions on the plane of the sky at the three epochs 
should lie along a straight line. A quantitative estimate of 
this is obtained from the correlation coefficient, r, of each 
triplet of positions per spot. We have hence computed r 
for all spots and calculated the arithmetic mean. This is 
shown in Fig. ^ as a function of the distance, d. A max- 
imum is reached for d=lA kpc, but also larger distances 
cannot be excluded, whereas for d < 1 kpc the correlation 
coefficient drops significantly. 

An alternative approach is that of comparing the mean 
direction of the proper motions (after subtraction of the 
apparent motion) with the direction of the jet axis. To 
give equal weights to all spots we have performed a vec- 

3 The contribution of the parallax and source motion rela- 
tive to the Sun have been calculated with a program kindly 
provided by Prof. Tetsuo Sasao. Both a flat galactic rotation 
curve and the rotation curve of Brand & Blitz (.19933 have been 
considered, but the results obtained for a distance of 1.7 kpc 
were not significantly different. 
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Table 1. Parameters of the water maser spots in IRAS 20126+4104. 



Number 


. . / „ \ 

ARA (a) 


. m 

ADec (a) 


s v 


Vlsr 


7T7\ 

Vra 


7T7\ 


ZrZ TT\ 


ZrZ m 




(mas) 


(mas) 


(Jy) 


(km s _1 ) 


(km s _1 ) 


(km s _1 ) 


(km s _1 ) 


(km s _1 ) 


1 


0.0+0.3 


0.0+0.3 


24.30 


-15.0 


-75+10 


17+10 


-84+10 


54+10 


2 


-1.0+0.3 


0.5+0.3 


18.36 


-15.3 


— 


— 


— 


— 


3 


-2.3+0.3 


1.2+0.3 


10.94 


-16.2 


— 


— 


— 


— 


4 


-1.6+0.3 


0.8+0.3 


5.12 


-16.0 


— 


— 


— 


— 


5 


-24.5+0.3 


-51.5+0.3 


4.50 


-2.3 


-39+11 


-12+10 


-48+11 


25+10 


6 


-21.5+0.3 


-65.8+0.3 


5.58 


-3.3 


-57+11 


-23+11 


-66+11 


14+11 


7 


-21.1+0.3 


-66.6+0.3 


0.92 


-5.1 


-48+12 


-23+11 


-58+12 


14+11 


8 


-22.9+0.3 


-63.0+0.3 


0.47 


-5.8 


-42+11 


-19+10 


-52+11 


19+10 


9 


-28.0+0.3 


-60.4+0.3 


5.77 


-8.5 


-46+11 


-20+11 


-55+11 


17+11 


10 


-25.8+0.3 


-61.6+0.4 


2.97 


-9.5 


-43+12 


-33+11 


-53+12 


4+11 


11 


-26.0+0.3 


-60.9+0.3 


2.08 


-8.6 


-43+11 


-21+11 


-52+11 


16+11 


12 


-9.2+0.3 


0.1+0.3 


0.29 


-14.8 


-67+11 


19+10 


-76+11 


56+10 


13 


12.3+0.3 


-83.8+0.3 


0.24 


-5.6 


-33+11 


-31+11 


-42+11 


7+11 


14 


11.7+0.3 


-84.9+0.3 


0.08 


-3.7 


-35+12 


-20+11 


-44+12 


17+11 


15 


14.8+0.4 


24.4+0.4 


0.49 


-8.0 


-19+12 


-19+11 


-28+12 


18+11 


16 


-141.9+0.3 


73.0+0.3 


0.85 


-24.9 


-76+11 


9+10 


-85+11 


46+10 


17 


-142.6+0.3 


73.2+0.3 


0.09 


-24.8 










18 


-17.8+0.3 


-67.0+0.3 


0.17 


-9.5 


-61+11 


-25+10 


-70+11 


12+10 


19 


-19.0+0.3 


-66.4+0.3 


0.14 


-7.9 


-49+11 


-28+10 


-59+11 


9+10 


20 


2.1+0.3 


-0.7+0.3 


0.13 


-14.2 


-66+11 


46+10 


-76+11 


83+10 


21 


6.4+0.3 


29.6+0.3 


0.12 


-8.0 










22 


29.9+0.3 


-77.1+0.3 


0.16 


-3.8 










23 


28.9+0.3 


-78.2+0.3 


0.13 


-3.8 










24 


21.3+0.3 


-87.2+0.3 


0.06 


2.7 










25 


28.2+0.3 


-78.0+0.3 


0.10 


-3.8 










26 


482.7+0.3 


-228.3+0.3 


1.06 


-1.9 


7+11 


-48+10 


-2+11 


-11+10 



(Q) offset with respect to RA(J2000)=20 h 14 m 26=0253, Dec(J2000)=41° 13'32'.'666 
absolute proper motion measured with phase referencing 

absolute intrinsic proper motion after correction for parallax, motion of the Sun with respect to the LSR, and galactic 
rotation 



tor average of the unity vectors of the proper motions in 
columns 6 and 7 of Table ^ weighted by the correspond- 
ing errors. The result is shown in Fig. ^ as a function of 
the distance. The grey area in the same figure marks the 
range of plausible values of the jet position angle (P.A.), 
obtained from tracers such as the SiO(2-l) and H 2 lines 
(Cesaroni et al lH^ . and the 2.2 fim (Edris et al. 00051 
hereafter EFCE) and 3.6 cm continuum emission (Hofner 
et al. 119991 Hofner pers. comm.). Clearly any distance 
between 1.2 and 2.7 kpc is possible. 

We conclude that with the present data only loose con- 
straints can be set on the distance to IRAS 20126+4104, 
indicating that the source cannot be closer than 1.2 kpc. In 
the following we will assume the usual estimate of 1.7 kpc 
for consistency with previous studies of the source. 

3.3. Distribution and velocities of the H2O maser spots 

As previously mentioned, the main scope of our study was 
to investigate the geometry and kinematics of the jet in 
IRAS 20126+4104 using the H2O maser spots as test par- 
ticles of the outflowing molecular gas. Figure^ shows all 
the H 2 maser spots detected in the region by us, MCR, 



1 1 r 




0.7 U , , , , 1 , , ^ 1 , , , , L 

12 3 



d (kpc) 

Fig. 1. Plot of mean correlation coefficient of the maser 
spots (dashed curve; see text) and mean direction of the 
normalised proper motions (solid curve) versus source dis- 
tance. The vertical dotted and solid lines mark respec- 
tively d=lA and 1.7 kpc. The grey area denotes the range 
of possible jet position angles. 



and EFCE. Also shown are the OH maser emission peaks, 
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which according to EFCE mark the plane of the circum- 
stellar disk. We overlay the masers on a contour map of the 
3.6 cm continuum emission imaged by Hofner et al. (jl999), 
corresponding to the inner ionised part of the jet/outflow 
in IRAS 20126+4104. All of these tracers agree very well 
with the pattern of the conical jet model that will be dis- 
cussed in Sect. 01 and is already presented here for the sake 
of comparison. As expected, the vertex of the cone - i.e. 
the putative position of the YSO powering the jet - falls 
very close to the peak of the 3.6 cm continuum emission 
and in between the two OH maser peaks denoting the disk. 

All this strongly favours the scenario already outlined 
by MCR, according to which the H2O masers are located 
at the interface between a conical jet and the surrounding 
quiescent material, and are flowing along the surface of 
such a jet with velocities directed outward from the cen- 
tral YSO. Such an interpretation is fully confirmed by the 
proper motions of the H2O maser spots, shown in Fig. [2b: 
clearly, the spot velocities are diverging from the YSO. We 
stress again that in order to obtain the intrinsic proper 
motions of the spots (columns 8 and 9 of Tabled the ab- 
solute proper motions (columns 6 and 7) have been cor- 
rected for the annual parallax and motion of the source 
relative to the Sun, assuming a distance of 1.7 kpc. 

Combining the proper motions with the LSR velocity 
of each spot after subtraction of the systemic velocity of 
-3.5 km s , one can obtain the corresponding speed of 
each spot. In Fig.|3we plot this (filled circles) as a function 
of the distance r (measured on the plane of the sky) from 
the peak of the 3.6 cm continuum emission, which is a 
good guess for the location of the YSO. Although with a 
large scatter, the speed tends to increase with r. A similar 
trend is seen for the l.o.s. components of the spot velocities 
- which represent a lower limit to the full speed. These are 
shown in the same figure (empty circles) also for the spots 
that do not have measured proper motions. We conclude 
that the jet is accelerated within at least 700 au from the 
YSO. 

For a more quantitative description of the maser dis- 
tribution and velocities, in the next section we elaborate 
a simple model which represents an improvement on the 
model developed by MCR. 



r (AU) 



300 



400 



500 



600 



700 



100 



50 



200 



250 



300 



350 



400 



(mas) 



Fig. 3. Plot of speed versus separation (measured on the 
plane of the sky) from the 3.6 cm continuum peak for all 
H2O maser spots of Fig. [2b with measured proper motions 
(filled circles). Also shown is a plot of the l.o.s. component 
of the velocity for all spots detected in this study (empty 
circles) . Note the tendency of the speed to increase with 
distance from the YSO. 



on the vertex, with z along the l.o.s. (the observer lies 
at z — —00) and x coincident with the projection of the 
jet axis on the plane of the sky. For further details and 
a sketch of the model we refer to Sect. 4.1 and Fig. 4 of 
MCR. 

Under the previous assumptions one can express the 
components of the velocity along the axes as 



dR 

— v 

dR y 

dv 

dR * 



(2) 
(3) 
(4) 



4. A model for the H2O masers in 
IRAS 20126+4104 

The scope of the jet model presented here is to prove that 
the scenario outlined in the previous section for the origin 
of the H2O maser emission in IRAS 20126+4104 is correct. 
Unlike MCR's model, ours takes into account that the 
spot speed increases with distance from the YSO and fits 
also the velocity components in the plane of the sky. We 
assume that the maser spots lie on the surface of a cone 
with semi-opening angle inclination angle with respect 
to the l.o.s. if), and vertex coincident with the YSO at 
aoA- The spot velocities are directed radially outward 
from the vertex and are proportional to the distance, R, 
from it: v(R) = R. The coordinate system is centred 



where ^ is constant and 



x sin ip cos ip (1 + cot 2 



± tan 9 J x 2 + y 2 (cos 2 if) - cot 2 9 sin if)) 



(5) 



x ( cot 2 9 sin 2 



if) — cos 2 if)) 



Note that Eq. (2) of MCR contains an unfortunate typing 
error (tan 2 9 must be replaced by cot 2 9) which makes it 
different from our Eq. JSJ. The correct equation was used 
in the calculations of MCR. 

The input parameters of the model are: the P.A. of 
the jet, 9, if), and the coordinates of the cone vertex, 
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Fig. 2. a) Map of the H2O maser spots detected towards IRAS 20126+4104 by us (small circles), MCR (squares), and 
EFCE (triangles) overlayed on a contour map of the 3.6 cm continuum emission (Hofner et al. I1999|l . Also shown are 
the OH maser emission peaks observed by EFCE (big circles). The errors on the absolute positions of the H 2 maser 
spots detected by us are given in Tablc^ while those of MCR and EFCE are respectively 30 mas and 15 mas; for the 
OH masers the error is 25 mas. The colour denotes the LSR velocity of each spot according to the colour scale in the 
bottom panel. The solid lines indicate the conical jet which represents the best fit to the H2O maser spots obtained 
with the model discussed in Sect. b) Enlargement of the central region illustrating the locations and absolute 
proper motions (corrected for parallax, solar motion with respect to the LSR, and galactic rotation) of the H2O 
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a and <5 • We searched for the values of these parameters 
minimising the expression 



26 
i=l 



y v 



+ (v®-VU) (6) 



where i indicates the spot number in Tableland symbols 
v and V denote respectively the model and observed ve- 
locities; for spots without measured proper motions only 
the term involving the z components is considered in the 
summation. 

The best fit is obtained for P.A.=123°, ip=96°, 9=17°, 
a (J2000)=20 h 14 m 26?0410, <5 (J2000)=41 o 13'32"536, and 
||=0.255 km s" 1 mas -1 =0.150 kms _1 au -1 . Figure H 
compares the measured values of the proper motions and 
LSR velocities (filled circles) with those obtained from 
the model fit (empty circles). Clearly the agreement is 
very good, proving that our model provides a satisfac- 
tory interpretation for the origin of the water masers in 
IRAS 20126+4104. 

A comparison between our best fit parameters and 
those of MCR reveals significant differences, the most 
striking being the different expansion velocities: 23 km s _1 
in MCR and 34-112 km s _1 in our case (see Fig.^J). This 
is due to the fact that the spots are moving almost in 
the plane of the sky, so that it is very difficult to correct 
for the inclination angle without proper motion informa- 
tion. As a matter of fact, ip is larger in our model than 
in MCR. Moreover, our value (ip=96°) suggests that the 
north-western lobe of the jet is pointing away from us, 
whereas in MCR (t/>=59°) the opposite occurs. Our re- 
sult is consistent with recent findings: on the one hand, 
the jet is very close to the plane of the sky on scales as 
large as 10", or 0.08 pc (Cesaroni et al.|T999); on the other 
hand, the jet axis is known to undergo precession, with the 
NW lobe receding from the observer at the present time 
(Cesaroni et al. I2005f) . It is hence not surprising that on 
scales as small as a few 100 au the jet has already crossed 
the plane of the sky, as indicated by the value ■0=96°. 

As for the opening angle of the cone, MCR's estimate 
(0=29°) is twice as much as ours (17°). This is partly re- 
lated to the fact that also the vertex of the cone has moved 
by ~84 mas to the SE with respect to MCR's estimate. As 
discussed in Sect. l3.3l the new position is closer to the disk 
plane, identified by the OH maser emission (see Fig. 

In conclusion, we believe that the estimates of the H2O 
maser jet obtained in the present study represent an im- 
provement on those of MCR and are hence to be preferred. 



5. Variability of the H2O masers 

It is interesting to complement the high resolution 
study of the spatial distribution of the H2O masers in 
IRAS 20126+4104 with an analysis of their time variabil- 
ity. Water masers are known to be highly variable and 
IRAS 20126+4104 is no exception to this rule. Since 1987, 
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Fig. 4. Comparison between the observed velocities (filled 
circles) and those obtained from our model fit (empty cir- 
cles). The top and middle panel show respectively the 
proper motions in the RA and Dec directions, while 
the bottom panel plots the LSR velocity. Note that 
data points missing in the top and middle panels corre- 
spond to spots without measured proper motions. The 
dashed line marks the systemic velocity (-3.5 km s _1 ) of 
IRAS 20126+4104. The number on the x-axis indicates the 
maser spot according to the notation in Table I5~3l 



the Medicina 32-m antenna 4 has been used as a single-dish 
telescope to monitor a number of galactic H2O masers 
with typical sampling intervals of 2-3 months. For de- 
tails on the observations and data analysis we refer to 
Valdettaro et al. ifgUoljl . IRAS 20126+4104 belongs to this 
sample, so that we could retrieve the H2O spectra from 
the Medicina database and analyse them with the method 

EE 



developed by Valdettaro et al. I)2()()2|l . This is shown in 
Fig. which represents an image of the maser intensity 
as a function of time and LSR velocity. Also shown are 
the points corresponding to the peak velocities of the H 2 
lines detected at different epochs, including those obtained 
by us with global VLBI, by MCR with the VLBA, and by 
EFCE with MERLIN. 



The Medicina VLBI radiotelescope is operated by the 
INAF-Istituto di Radioastronomia. 
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VLBA VLBI MERLIN 




t (days) 

Fig. 5. Variability plot of the H2O maser emission in IRAS 20126+4104. The grey scale map represents the intensity 
measured with the Medicina antenna as a function of time and LSR velocity. The value t=0 corresponds to September 
12, 1995. The white squares mark the positions of the the peak velocities of the lines detected in the Medicina spectra, 
while the error bars correspond to the line full width at half maximum (FWHM). Both the peak velocity and line 
width have been obtained with a Gaussian fit. The black squares indicate the peak velocities of the spots observed in 
the VLBA, VLBI, and MERLIN observations of MCR, this study, and EFCE respectively. Only the first and third 
VLBI epochs are shown here to preserve readability of the figure. The horizontal solid line marks the systemic velocity 
(-3.5 km s _1 ) of IRAS 20126+4104. The dashed lines outline the linear variation of velocity versus time observed in 
three spectral features, labeled with numbers (1) to (3). 



The most striking characteristic is the presence of three 
clear velocity trends, outlined by the dashed lines in Fig. [5] 
in all cases, the velocity decreases linearly with time, with 
slopes dVLsR/di - -2.1,-2.7, and -2.1 km s^yr" 1 for 
features (1), (2), and (3) respectively. In the following we 
discuss these features in detail. 



5.1. Features (2) and (3) 

For features (2) and (3) one may easily derive the value 
of the deceleration corrected for the inclination angle. In 
fact, we have measured the proper motion of the corre- 
sponding spots: using the notation of Tabic ^ these are 
spots 1, 12, and 20 for feature (2) and 16 for feature (3). 
The inclination angle is given by the ratio between the 
l.o. s. velocity and the projection of the velocity on the 
plane of the sky. Applying this correction, one obtains 



dV/dt ~ -18 km s^yr -1 fOT ( 2 ) and -9.7 km s^yr -1 
for (3). 

In Fig. the distribution of the spots corresponding 
to feature (2) is shown, including our measurements and 
those of MCR and EFCE. Notwithstanding the large po- 
sitional errors on the VLBA and MERLIN measurements, 
the general trend follows the direction of the H2O maser 
jet depicted in our model. This is proved by the empty cir- 
cles, which mark the positions of the spots with measured 
proper motions (filled circles with arrows) extrapolated at 
the two epochs of the VLBA and MERLIN observations. 
Such an extrapolation takes into account the acceleration 
in the plane of the sky estimated above. Clearly, the lo- 
cation of the empty circles is consistent within the errors 
with that of the VLBA (squares) and MERLIN (triangles) 
spots. 

It is worth noting that the error bars in Fig.^represent 
the uncertainty on the absolute position at each epoch. 
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Fig. 6. Map of the H2O maser spots belonging to fea- 
ture (2) in Fig. |3] The (0,0) position corresponds to 
RA(J2000)=20 h 14 m 26 i !0253, Dec(J2000)=41°i3'32"666. 
Black circles indicate our measurements, squares those by 
MCR, and triangles those by EFCE. The error bars in- 
dicate the Id uncertainty in the astrometry of the MCR 
and EFCE measurements (the error on our data is less 
than the size of the symbols). Note that the error in the 
relative positions of the spots in each data set is of order 
1 mas. The arrows denote the proper motions measured 
by us, while the empty circles mark the positions of the 
same spots extrapolated to the dates of MCR (November 
21, 1997) and EFCE's (March 2002) observations, on the 
basis of the current spot position, proper motion, and ac- 
celeration (see text). 

The relative position of the spots observed at a given epoch 
is instead very accurate, of order ~1 mas. Therefore, it 
is clear that in the first epoch (squares) the dispersion 
between the spots was much greater than at the time of 
our VLBI observations. One may speculate that this an 
indication of the fact that the shocks associated with the 
masers are converging. 

5.2. Feature (1) 

Feature (1) has not been detected in our VLBI observa- 
tions, so that no proper motion measurement is available. 
However, the VLBA and MERLIN measurements of MCR 
and EFCE can be used to obtain an estimate of the aver- 
age proper motion. Notwithstanding the large uncertainty 
on the absolute positions of the spots, the time interval 
between the two measurements is large enough to guar- 
antee an angular separation between the positions at the 
two epochs much larger than the errors. This can be seen 
in Fig. [3 where also a sketch of the conical model pre- 
sented in Sect. 21 is shown for the sake of comparison. The 
mean velocity in the plane of the sky is obtained from 
the separation between the two positions (169±45 mas 
or 287±77 au) divided by the corresponding time inter- 
val (~1620 days), and turns out to be 306±42 km s _1 . 
We also stress that the direction of the proper motion is 





\f 1 1 1 1 1 


1 1 1 1 1 1 1 1 
1 — i 


» — 1 ; 


f 


— i 


1 1 




(1) ; 



150 100 50 

Aa (mas) 



Fig. 7. Map of the H2O maser spots belonging to fea- 
ture (1) in Fig. [SJ The (0,0) position corresponds to 
RA(J2000)=20 h 14 m 26?0253, Dec(J2000)=41°13'32"666. 
Squares indicate measurements by MCR and triangles 
those by EFCE. The error bars indicate the la uncer- 
tainty on the astrometry. The error on the relative posi- 
tion of the spots in each data set is of order 1 mas. The 
starred polygon marks the position of the YSO according 
to our conical model fit, while the solid lines denote the 
projection of the cone on the plane of the sky. 

perfectly consistent with expansion along the cone as ex- 
pected in our model: this has the twofold consequence of 
confirming our model and proving that all spectral fea- 
tures falling on the same dashed line in Fig. [5] do arise 
from a given spot (or group of spots) moving along the 
jet. 

One may compute the mean velocity of feature (1) 
along the l.o.s. from the arithmetic mean between the l.o.s. 
velocities at the two epochs (~11.2 km s _1 ). Finally, from 
the ratio between this and the proper motion, it is possible 
to obtain the inclination of the velocity vector with respect 
to the l.o.s. and thus correct the l.o.s. values of the mean 
speed (11.4 km s _1 ) and acceleration (-2.1 km s _1 yr -1 ): 
the corrected values are 306 km s _1 and -56 km s" 1 yr~ x . 

5.3. Origin of the deceleration 

From the previous findings, we conclude that H2O masers 
in IRAS 20126+4104 are undergoing a strong deceleration 
ranging from -50 to -10 km s _1 yr _1 . Although evidence 
of deceleration from single-dish monitoring of H2O masers 
has already been reported by other authors (see Brand et 
al.|5J)03), to our knowledge this is the first time that such 
a finding is confirmed by high-angular resolution observa- 
tions. 

This existence of deceleration may look contradictory 
with the positive velocity trend shown in Fig [21 However, 
dV/dR > does not necessarily imply dV/dt > 0, be- 
cause the V(R) provided by our VLBI measurements cor- 
responds to an instantaneous picture of the spots velocity 
field, whereas V(t) obtained from the variability analysis 
describes the evolution of each single spot with time. 

In order to reconcile dV/dR > with dV/dt < 0, 
we propose the following scenario. Water masers are very 
likely generated in the post-shock region of J-type shocks 
(Elitzur et al. 1989 ). In our model, such shocks occur along 
the conical surface representing the interaction between 
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Fig. 8. Peak velocity of feature (1) and (2) (see Fig. [SJ 
as a function of time. Error bars denote line FWHMs. 
The dashed lines are the same as in Fig. [S] and represent 
linear fits to the data, while dotted and solid lines are fits 
of the type V oc (t — to) -1 ^ 2 and V oc exp[-(t — t )/r], 
respectively. 



the jet and the surrounding cloud. As predicted by the- 
oretical models (e.g. Ferreira & Casse 2004, Shu ct al. 
2000}, the jet should be accelerated close to the YSO. 
Consequently, the further away (from the YSO) the jet 
impinges against the surrounding material, the larger will 
be the velocity of the corresponding shock. At this time 
a maser spot "is switched on" , with the same velocity as 
the shock. This explains why dV/dR > 0. After that, the 
shock is slowed down as it proceeds through the dense ma- 
terial of the cloud, thus causing the observed deceleration, 
dV/dt < 0, of the spot velocity. 

It remains to be verified that the scenario depicted 
above can fit the observed trend of V(t) quantitatively. 
If the deceleration is due to braking of a shock, this 
should imply a power-law for V(t). In fact, for shocks 
conserving energy one has MV 2 =const. =>■ yiV^ 2 =const. 
+> V oc t -1 ' 3 , while shocks conserving momentum im- 
ply MV=const. and hence V oc t~ x / 2 . Another possibil- 
ity is that of an exponential decay of the shock velocity: 
naively, this may be described with a sort of "viscosity" 
in the medium ahead of the shock. With this in mind, 
we have fitted the velocity trends of features (1) and (2) 
with V oc (t — to)^ 1 / 2 and V oc exp[— (t — io)/V]: these are 
represented respectively by the dotted and solid curves in 
Fig. |H1 where we plot the l.o.s. velocities of features (1) 
and (2) as a function of time. Also shown are linear fits 
to the data (dashed lines). It is worth noting that in all 
cases the lifetime of the maser spots ranges between 5 
and 12 yr, a timescale comparable to our monitoring pe- 
riod. Clearly, it is difficult to rule out any of the models, 
although the exponential fit seems more satisfactory. A 
more sensitive and frequent monitoring of the H2O maser 
emission in IRAS 20126+4104 is necessary to discriminate 
between models. 



6. Summary and conclusions 

We have performed three-epoch VLBI observations of 
the H2O maser emission from the massive protostar 
IRAS 20126+4104. Phase referencing allows us to mea- 
sure the absolute proper motions of the maser spots and 
hence obtain a 3-D picture of the velocity field. These 
high angular resolution results have been complemented 
by an analysis of the variability of the H2O maser emis- 
sion through single-dish observations spread over 15 years. 
The main results obtained from this study and comparison 
with previous data by MCR and EFCE are the following: 

— A total of 26 spots have been detected in our VLBI 
observations, 17 of which at all three epochs thus al- 
lowing measurement of their proper motions. 

— An attempt to use proper motions to obtain an esti- 
mate of the distance to the source may only loosely 
constrain it between 1.2 and 2.7 kpc. 

— After correcting the measured proper motions for an- 
nual parallax and apparent motion of the source with 
respect to the Sun (assuming a distance of 1.7 kpc), 
we find that the maser spots are expanding from a 
common origin, consistent with the expected position 
of the protostar. The speed is increasing from 34 to 
112 km s _1 for increasing distance from the protostar. 

— The l.o.s. velocities and proper motions of the spots 
can be fitted with an improved version of the coni- 
cal jet model proposed by MCR in a previous VLBA 
study of the same maser source. The best fit is ob- 
tained for a jet semi-opening angle of 17°, an incli- 
nation of 96° with respect to the l.o.s. (the jet axis 
to the NE is pointing away from the observer), an 
expansion velocity gradient of 255 km s _1 yr _1 , and 
a cone vertex located at RA(J2000)=20 h 14 m 26!0410, 
Dec(J2000)=41°13'32'.'536. This is consistent with the 
expected position of the protostar. 

— The variability study demonstrates that some of the 
maser features undergo decelerations from ~ -50 to 
-10 km s _1 yr _1 , with lifetimes of the order 5-12 yr. 
We speculate that this may be due to braking of the 
shocks where the maser emission is originating, at the 
interface between an accelerated jet and the surround- 
ing molecular environment. 

We conclude that our study has shed light on the 3-D 
structure and kinematics of the jet in IRAS 20126+4104 
and outlined some intriguing characteristics of the time 
variability of the maser - and hence of the jet. A full un- 
derstanding of this latter aspect will require more sensitive 
observations with better angular and, possibly, temporal 
resolution. 
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